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Abstract 

The estimation of the solution of a system of two differential equations introduced 
by Norton et al. (1976) that is equivalent to the famous Gompertz growth law is 
performed by means of the recent adaptive scheme of Besangon and collaborators 
(2004). Results of computer simulations illustrate the robustness of the approach. 
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1 Gompertz growth functions 



Mathematically speaking, the Gompertz law refers to the class of functions 
having exponentially decreasing logarithmic derivatives. It has been intro- 
duced in 1825 by B. Gompertz in a paper on the law of human mortality (T). 
He noted that for people between 40 and 100 years "t/ie rate of mortality vari- 
able with age measures a force of death which grows each year by a fraction, 
always the same, of its present value". According to Winsor i^, the possible 
application of the Gompertz curve in biology was first spelled out in 1926 by 
Sewall Wright in a few remarkable statements: " In organisms, on the other 
hand, the damping off of growth depends more on internal changes in the cells 
themselves ... The average growth power as measured by the percentage rate 
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of increase tends to fall at a more or less uniform percentage rate, leading to 
asymmetrical types of s- shaped curves...". 

If the size z{t) of a growing structure evolves according to the equation (j^ 
z\n(^] , (1) 



we say that its growth is of Gompertz type. The evolution is continuous from 
a given initial stage to a plateau value S. In a Nature letter on the growth 
of tumours, Norton et al ^ wrote the Gompertz law as the system of the 
following two first order differential equations 

Zi = K,ZiZ2 (2) 

Z2 = -K2Z2, (3) 



where Z = {Zi, Z2) G M^, Ki > 0, Zi is the volume of the tumour at time t, 
and Z2 is a function entirely described by the second equation (3) that gives 
the difference in growth with respect to a pure exponential law. According 
to Norton, K2Z2 gives the fraction of the volume that doubles in size during 
the instant dt. Thus, Z2 that we call for obvious reasons the Gompertzian 
exponent function is of special interest and we would like to determine it with 
high accuracy even though we know neither the initial conditions for Zi and 
Z2 nor K2. Norton et al wrote the solution of the system in the following form 



Zi(t) = Zi(0)exp 
^2(t)=(^Un 



K2 
S 



Z2(0){l-exp(-ir2t)} 



Ziit) 



(4) 
(5) 



We will treat Zi and Z2 as states of a dynamical system that in our case is 
the evolution of a tumour. The fundamental concept of state of a system or 
process could have many different empirical meanings in biology and in our 
case the first state Zi is just the size of the tumor whereas Z2 is the deviation 
of the growth rate from the pure exponential growth. In general terms, a 
potentially useful tool in biology is the reconstruction of some specific states 
under conditions of limited information. For animal tumors, it is not trivial 
to know their state Zi at the initial moment and most often we do not know 
the instant of nucleation that can be determined only by extrapolation of the 
fitting to the analytic solutions of growth models, such as Eqs. (4) and (5). In 
this paper, the main goal is to show that an excellent alternative procedure 
for estimating the phenomenological quantities of the tumor growing process 



2 



in the frequent case in which we do not know the initial conditions and the 
parameter K2 is the recent adaptive scheme for state estimation proposed 
by Besangon and collaborators p|). In addition, what is generally measured, 
i.e., the output is a function of states that we denote by h{Z) and in the 
particular case of tumors one usually measures their volume. Then: 



y = h{Z) = Z, . (6) 

The mathematical formalism of the adaptive scheme that follows relies entirely 
on the Lie derivatives of the function h{Z) that are defined in the next section. 
By a Lie mapping, we are able to write the Gompertz-Norton system in Be- 
sangon's matrix form (system F below) that allows to write the corresponding 
adaptive algorithm (the dynamical system F and its explicit Gompertz form 
F G below). 



2 The adaptive scheme 



Taking into account the fact that rarely one can have a sensor on every state 
variable, and some form of reconstruction from the available measured out- 
put data is needed, an algorithm can be constructed using the mathematical 
model of the process to obtain an estimate, say X of the true state X. This 
estimate can then be used as a substitute for the unknown state X. Ever since 
the original work by Luenberger (0), the use of state 'observers' has proven 
useful in process monitoring and for many other tasks. The engineering con- 
cept of observer means an algorithm capable of giving a reasonable estimation 
of the unmeasured variables of a process using only the measurable output. 
Even more useful are the so-called adaptive schemes that mean observers that 
are able to provide an estimate of the state despite uncertainties in the pa- 
rameters. The so-called high gain techniques proved to be very efficient for 
state estimation, leading in control theory to the well-known concept of high 
gain observer (0). The gain is the amount of increase in error in the observer's 
structure. This amount is directly related to the velocity with which the ob- 
server recovers the unknown signal. The high-gain observer is an algorithm 
in which the amount of increase in error is constant and usually of high val- 
ues in order to achieve a fast recover of the unmeasurable states. In case of 
dynamical systems depending on unknown parameters, the design of the ob- 
server has to be modified appropriately in order that the state variables and 
parameters could be estimated. This leads to the so called adaptive observers, 
i.e., observers that can change in order to work better or provide more fit for 
a particular purpose. Recently, observers that do not depend on the initial 
conditions or the estimated parameters from the standpoint of asymptotic ex- 
ponentially fast convergence to zero of the errors have been built for many 
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systems. They are called globally convergent adaptive observers and have been 
obtained from a non trivial combination of a nonlinear high gain observer and 
a linear adaptive observer, see (8) and (j^)- In this work, we present an appli- 
cation of the high gain techniques in the context of state estimation whatever 
the unknown parameter is. 

The assumption on the considered class of systems are basically that if all of 
the parameters were known, some high-gain observer could be designed in a 
classical way, and that the system are "sufficiently excited" in a sense which is 
close to the usually required assumption on adaptive systems, that is, signals 
should be dynamically rich enough so that the unknown parameters can indeed 
be identified. In this particular case, the lack of persistent excitation of the 
system could impede the reconstruction of the parameters. However, the recent 
scheme of Besangon and collaborators (0) guarantees the accurate estimation 
of the states according to rigorous arguments in their paper. 

To make this mathematically precise we have to introduce first some termi- 
nology. Let us construct the jth time derivative of the output. This can be 
expressed using Lie differentiation of the function h by means of the vector 
field / given by the right hand sides of Norton's system. We will denote the jth 
Lie derivative of h with respect to / by Lf^ (h) {Z{t)). These Lie derivatives 
are defined inductively as functions of Z 



L/ih) iZ) = ^(Lr'ih) {Z))f{Z). 



Lf\h) {Z) = h{Z) 
d_ 
'dZ 

When the system is observable, i.e., from the knowledge of the output one can 
build the states of the system, the Lie map $ : Z — $(2') given by 



( HZ) \ _ 




{LfihiZ))) 


K^Z^Z2 



(7) 



is a diffeomorphism. For ^{Z) to be a diffeomorphism on a region Q, it is 
necessary and sufficient that the Jacobian d$(Z) be nonsingular on Q and 
that ^{Z) be one-to-one from Q to see 

Since ^{Z) is a diffeomorphism, one can write the global coordinate system 
defined by X = $(Z) in the following form 
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T : 



Xi — X2 

,2 



X2 = ^ - K2X2 

y = x. 



Following ^ , we assume that the T system can be written in the matrix form 
as follows 



X = AX + ^{x) + -^{x)e 

y = CX 



(8) 



where X G 



A 



1 




, ^(X) = (0, 



Xi 



y is the measured output, \E'(t) G 



j^nxp matrix of known functions and G is the vector of unknown 

parameters that should be estimated through the measurements of the output 
y. We are here in the particular case n = 2, p = 1, i.e., \E'(t) = [0, — X2]"^ and 
G = K2. In addition, the algorithm we develop is a particular case of 
that presented in (|5|), since we can not meddle in the system, in other words, 
there is no control input. As adaptive observer we use the system 



X = AX + ip{a{X)) + ^{a{X))a{e) + [p/C + TT^C^] {y - CX 



F = < 



e ■■ 
r 



CX 



p(A-/cc)r + p^(cT(x)) 



(9) 



where a(-) is a saturation function, F G 



Dnxp 



is the so-called gain vector, /C is a 



vector that makes A — ICC a stable matrix, A = diag[l, p 



"("■ 1)1 where 



p G IR+ is a constant to be chosen. The saturation function is a map whose 
image is bounded by chosen upper and lower limits, B and b, respectively. It 
is customary to introduce such functions of simple forms, e.g., we used 



as 



B s> B 
<; s b< s <B 
) s < b , 



(10) 



to avoid the over and/or underestimation and in this way to increase the 
chance of the quick convergence to the true value (fl^ . 

In (j^, it is proven that the dynamical system F is a global exponential adap- 
tive observer for the system F, i.e., for any initial conditions X(to), X(to), 
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0{to) and we E W, the errors X{t) - X{t) and e{t) - 9{t) tend to zero expo- 
nentially fast when t — > oo. Taking JC — [/Ci,/C2], the matrix A — ICC have the 
following eigenvalues 

Ai,2 = -1/2 /Ci ± 1/2 ^kJ^AIC^ . (11) 



Selecting IC2 = (1/4) /Ci^, we get equal eigenvalues Ai = A2 = — (1/2) JCi, and 
choosing /Ci > we turn A — K.C into a stable matrix. Thus, the explicit form 
of the observer system F is given by 



(^(^2)) 



a(X2)a(^^) + p + rir2) {Xi - X,) 



9 = p^T^iX, - X^) 

fi = p(~/Ciri + r2) 

r2 = -£^ + pa(X2) 



(12) 



Being global, this observer system does not depend on the initial conditions. 
Therefore, any initial conditions chosen at random from a set of physical val- 
ues will not affect the correct estimation; merely the convergence time could 
be longer or shorter. Thus, in practice, it is useful to start with initial condi- 
tions that are close to the real phenomenological initial conditions in a given 
framework. 

Finally, to recover the original states, we use the inverse transformation $~^(X), 
which is given by: 



Xi 



(13) 



With the aim of better illustrating the adaptive scheme proposed here, we 
present numerical simulations. We use the following values of the parameters: 
Ki = 1, K2 = 0.5, p = 100 and /Ci = 1. In Figs. (1) and (2), the solid lines 
represent the evolution of the true states and the dotted lines stand for the 
evolution of the estimates, respectively. We mention that short convergence 
time is what really matters in order to have efficient numerical simulations. 
This can be accomplished by starting with arbitrary initial conditions that 
are guessed to be close to the real initial ones as we already commented. If 
one is interested in the evolution of the iterative scheme, this can be readily 
glimpsed from the difference Z^ — Zi between the curves in the figures. 
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To this end, we would like to illustrate the robustness of the present adaptive 
scheme. Figs. (3) and (4) show what happens when an impulsive type pertur- 
bation (i.e., of high value acting in a very short span of time) is added to the 
output signal h{Z) = Zi which is fed to the observer at t = 4 (arbitrary units). 
This perturbation is a strong error in the output y that can be attributed to 
an unexpected fault in a measurement device or even human errors and does 
not belong to the natural evolution of the tumour. As can be seen from the 
graphics, the adaptive scheme has the ability to recover the "true" signal im- 
mediately after the perturbation disappears. In general, this robustness is due 
to the fact that the scheme is designed in the closed-loop way and additionally 
not the full range of the parameters need to be known. 



3 Conclusion 

The robust adaptive scheme we used here for the interesting case of Gompertz 
growth functions is a version of that due to Besangon et al. The results of 
this work indicate that this scheme is very efficient in obtaining the Gompertz 
functions without knowing both initial conditions and parameter K2- The 
method may be useful in more general frameworks for models of self-limited 
growth such as in the construction of a specific growth curve in biology, or 
as a managerial tool in livestock enterprizes, as well as in the detailed under- 
standing of the growth of tumors. We also notice that the reconstruction of 
the unknown states by this method allows the possibility to obtain important 
missing parameters by standard fitting procedures. 
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Fig. 1. Numerical simulation for the first state: (a) the solid line represents the 
time evolution of the true states Z\ and the dotted line represents the estimate Z\. 
Plot (b) is a detail of (a) to appreciate the variation of Z\ in the beginning. 




time time 

Fig. 2. Numerical simulation for the second state: (a) the solid line represents the 
true state Z2 and the dotted line represents the estimate Z2. Plot (b) is a detail of 
(a) in order to appreciate better the variation of Z2 in the beginning. 
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Fig. 4. Behaviour of Z2, Z2 under the same perturbation at t = 4 (arbitrary units). 
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